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Abstract 


In this study, a two-fluid (TF) model is developed for two-phase flows in proton exchange membrane fuel cells (PEMFCs). The drag force and lift 
force between gas and liquid phase are considered in N-S equations. In addition, a simplified model is introduced to obtain the liquid water droplet 
detachment diameter on the gas diffusion layer (GDL)/channel interface which involves the properties of the GDL/channel interface (contact angle 
and surface tension). The TF model and the simplified model for the prediction of water droplet detachment diameter on GDL/channel interface 
are validated by the comparison between the experimental data and the model results, respectively. The effect of the properties of GDL/channel 
interface (contact angle and surface tension) on two-phase behavior in PEMFCs is investigated, The results show that a high contact angle and a 
low surface tension are advantageous for liquid water removal in the gas channel and the GDL even though a low surface tension will lead to a low 


capillary force in the GDL. 
© 2006 Elsevier B.V. All rights reserved. 
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1. Introduction 


Proton exchange membrane fuel cells (PEMFCs) are 
regarded as the most promising energy conversion systems for 
future automobiles and stationary applications nowadays. Water 
management is the key issue in PEMFCs, and is a significant 
technical challenge. On the one hand, sufficient water is needed 
in the membrane to maintain sufficiently high proton conductiv- 
ity, but on the other hand, excess liquid water in the electrode can 
cause ‘flooding’, and hinders the transport of the reactant from 
the gas channel to the catalyst layer. Liquid water transport in 
PEMEFCs occurs as follows: (1) water is produced in cathode cat- 
alyst layer, and liquid water transports within the gas diffusion 
layer (GDL) by capillary-driven flow. (2) Liquid water droplets 
appear on the GDL/gas channel interface and are removed by the 
gas shearing function [1—4]. (3) Liquid water travels in the gas 
channel. Wang and co-workers [3] have observed the emergence, 
growth and detachment of liquid water droplet on the GDL/gas 
channel interface. Recently, they [5] have measured the size of 
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the droplet on the GDL/gas channel interface, and find that liquid 
water can transport through the gas channel without interaction 
with channel wall at high gas velocity. Modeling liquid water 
behavior in PEMFCs will help to optimize the water manage- 
ment and the design of PEMFCs. At present, much modeling 
work has been done [6-12]. He et al. [6] developed a two-phase 
model of water transport in GDL for a two-dimensional cross- 
section of PEMFCs. They introduced a transport equation for 
liquid water in GDL which involved an advection term (due to 
bulk flow) and a diffusion term (due to capillary-driven flow). In 
Siegel’s work [7], he assumed that the liquid phase had the same 
velocity with the gas phase, so the work is independent on the 
transport of liquid water in gas channel. He also considered the 
mass transfer in the catalyst layer in detail, which includes mass 
transfer among liquid phase, gas phase and polymer phase. Liu 
and co-workers [8—10] have employed the mixture model for 
two-phase behavior in PEMFCs, which can describe the liquid 
water behavior in GDL and gas channel with the advantages of 
simplicity and low cost of the calculation. Wang et al. [11] also 
introduced the mixture model to describe two-phase behavior 
in PEMFCs, and they predicted that the liquid water saturation 
within the cathode will reach 6.3% at 1.4 A cm? for dry inlet 
air. In the continuation of this work, Wang and co-workers [12] 
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Nomenclature 

C concentration (mol cm~?) 

CL lift coefficient 

CLa constant for Cz (0.2) 

de mean pore diameter 

D diffusion coefficient (cm? s7!) 
Dp characteristic droplet size (cm) 
F Faraday’s constant 

h water droplet height (cm) 

J reaction rate 

K constant for droplet size 

M molar mass (kg mol`!) 

M interface force of liquid (N) 
Mg interface force of gas (N) 

n number of electrons 

P pressure (Pa) 

q switch function 

R universal gas constant (J mol`! K7!) 
Re Reynold number 

sS water saturation 

T temperature (K) 

U velocity (cm/s) 

Veell cell voltage (V) 

Voc open circuit voltage (V) 

x mass fraction 

y molar fraction 


Greek letters 


E volume fraction 

(0 potential (V) 

À polymer water content H20/SO3 
u viscosity (Pa s) 

é equilibrium contact angle on diffuser 
p density (kg cm7?) 

o surface tension (N cm7!) 
Subscripts 

a anode 

c cathode 

C about capillary 

d diffusion layer 

Darcy Darcy pressure loss 

e about electron 

g gas phase 

H2 hydrogen 

i note for species 

j note for species 

l liquid phase 

m polymer phase 

mix gas mixture 

O2 oxygen 

R relative 

w water 

wd water dissolved in polymer 


wp water production 
wv water vapor 
Superscripts 

D drag force 

eff effective 

L lift force 

p polymer 

ref reference 

sat water saturated 


simulated the function of a microporous layer between cath- 
ode GDL and catalyst layer with the mixture model. Meng and 
Wang [13] have developed an experimental submodel to con- 
sider liquid coverage at the GDL/gas channel interface. In the 
present work, the size of liquid water droplet detached from the 
GDL/gas channel interface is considered, which has effect on 
the liquid removal in the channel and involves the properties of 
GDL/gas channel interface such as contact angle, surface ten- 
sion, average pore size of GDL and so on. A two-fluid (TF) 
model is employed to model liquid water droplet transport in 
the gas channel, the effect of GDL/gas channel interface prop- 
erties on the two-phase behavior in the gas channel and GDL is 
investigated. 


2. Model development 
2.1. Model assumptions 


The calculated regions consist of conventional channels, gas 
diffusion electrodes, catalyst layers and a membrane, which are 
shown in Fig. 1. Because PEMFC is operated in temperature 
below 100°C and above normal pressure, so it is assumed that 
water is generated in cathode catalyst layer as liquid. When 
water vapor concentration in GDL or gas channel is less than 
the saturated concentration, liquid water evaporates to vapor, 
otherwise, water retains a liquid state. The liquid water transfers 
within the GDL by capillary force, then the liquid water droplet 
emerges and grows on the GDL/gas channel interface. When 
the drag force acting on the droplet equal to the adhering force, 
the droplet obtains the maximum droplet size, then it begins to 
move on the GDL/gas channel interface or goes into the gas 
channel. The liquid water is treated as a continuous phase in the 
form of droplet. To hold up this approach, necessary assumption 
should be made, i.e. when the droplet detaches from GDL surface 
by gas shear flow, the effect of gas channel wall’s properties 
(i.e. hydrophobicity) are not considered. However, it is more 
reasonable for higher gas velocity or small droplet size, and 
Wang and co-workers [5] suggested that for a gas flow channel 
with the height of 0.5 mm, when the droplet diameter exceeds 
0.2 mm, the droplet may interact with gas flow channel walls. 
In addition, isothermal condition is also assumed to PEMFCs 
and gas mixture both in cathode and anode are treated as perfect 
gas. 
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Fig. 1. Calculated domain of PEMFCs. 


2.2. Governing equations 


A two-dimensional two-fluid model is used to describe two- 
phase flow behavior in a PEMFC gas channel. It considers the 
flow of gas and liquid phase, respectively. Terms accounting 
for the interacting force are added into the N-S equations. In the 
channel, the pressure of the two phases is assumed to be equal. In 
porous electrodes, the difference of the pressure depends on the 
capillary force, which is calculated from the Leverett J-function. 

Eqs. (1)-(12) in Table 1 contain the governing equations for 
two-dimensional two-fluid model of PEMFCs, which include 
conservation of mass, momentum, ionic charge and species. Eqs. 
(1) and (10) describe the conservation of mass of gas in cath- 
ode and anode, respectively. The mass conservation of liquid 
phase is described by Eqs. (5) and (7) for liquid in cathode gas 
channel and GDL, respectively. Eqs. (2), (6) and (11) describe 
the momentum transfer for gas phase in cathode, liquid phase 


Table 1 
Governing equations 


in cathode and gas in anode, respectively. The drag force and 
lift force between gas phase and liquid phase are considered 
in gas channel in Eqs. (2) and (6). Eqs. (3) and (12) describe 
the conservation of water vapor of cathode gas and anode gas, 
respectively. And the conservation of oxygen is given in Eq. 
(4). The conservation of ionic charge and water within polymer 
phase are expressed in Eqs. (8) and (9), respectively. Source 
terms reflect the mass exchange among liquid phase, gas phase 
and polymer phase as well as the consumption of reactant by 
electrochemical reaction are given in Table 2. Parameters used 
in the equations are listed in Table 3. 


2.3. The simplified model for calculating the droplet 
detachment diameter 


In PEMFCs, liquid water is produced in cathode catalyst layer 
and moves out of the gas diffuser by capillary-driven, and then 


Phase Conservation equation General form 

Cathode gas Mass V(e(l — 5)pgUg) = rw + So, — Swadwd (1) 
Cathode gas Momentum Vel — s)PgUgUg] = —e(1 — s)V pg + VIC. — S)UgV Ug] + EPg8g + SDarcy,g + (2) 
Cathode gas Species Vell — s)PgUgxo,] = DË V yo, + So, (3) 


Vie — 5)ogU exw] = DÉ V ywy + rw + qwa Swa 


(4) 


Cathode liquid (channel) Mass V(eso U1) = —rw (5) 
Cathode liquid (channel) Momentum V(espU1U1) = —esV py + V(esmV U1) + Sparey,1 + MP + M} (6) 
k 
Cathode liquid (electrode) Mass V(espU1) = V ( ESPI : V(Pg po) = —Tw + Swp — (1 — qwa)Swa (7) 
H 
J, anode catalyst layer 
Polymer phase Current V(KmV¢m) = 4 0 membrane (8) 
Je cathode catalyst layer 
A Swa 
2.5 e ue 
Polymer phase water Mass V(—DwVCw +a g tC,U)= (1 — qwa) M, catalyst layer ð 
-i 0 membrane 
Anode gas Mass V(epU) = Su, — qwa Swa (10) 
Anode gas Momentum V(epUU) = —eVp + V (eu VU) + £pg + Sdarey (11) 
Anode gas Species VlepUxwv] = V(oDẸ y, Vw) — qwa Swa (12) 
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Table 2 
Source terms in governing equations 
Source terms (zero in other region) Defining equation 
. Mo, 
Oxygen reaction rate (cathode catalyst layer) So, = —Je TF (13) 
. Mu, 
Hydrogen reaction rate (anode catalyst layer) Su, = —Ja OF (14) 
el—s 7 . 
Mass transfer rate between gas and liquid [14] rw = Mw [k ( = OwP peg + keesCw(Yw DP — p% —q) (15) 
Mass transfer rate between gas and polymer Swa hml o$, — py) (16) 
Darcy pressure drop of gas in cathode and anode Sparcy,g —& K 5 Ug (7) 
p(l -=s 
Darcy pressure drop of liquid in cathode SpDarcy,1 65 Ui (18) 
ps 
Liquid generation rate in cathode catalyst layer Swp = T Mw (19) 
Table 3 
Parameters in governing equations 
Parameters Defining equation 
Reaction rate in cathode Je = (1 — s)Ayittt Com [Fg ei (bm + Veet — Voc) re bm + Veen — Voe) (20) 
action rate in catho e = (1 —s)Ayi X| m + Vee oc ex] m + Vee oc 
0,c cg P RT 1l P RT 1l 
Reaction rate in anode J= Ai Crom e z (Pe — bm) = (Qe — bm) (21) 
action rate in 4 a Ay X e m ex] e k 
oa (e P (RT P\ RT 
1 1 
Polymer phase conductivity Km = (0.005139 — 0.00326) exp [1268 (= = z)| (A > 1) (22) 
Water content in polymer ph a (23) 
ater content in polymer phase a= 
poty P PET 
Water index in polymer à = 0.043 + 17.81a — 39.85a? + 36.0a7 24) 
EP Xwp Mmix 
Water activity = (25) 
Py M 
bs(p¥t — yw p + le — 30 
Liquid/vapor switch function q=0.5+4 Gin —Wwp tle ) (26) 
2(p\" — Ywp + le — 30) 
bs(p — ph + le — 30 
Gas/polymer switch function qaw = 9.5 4 ss ae 2 4 ) (27) 
2(pw — Pw + le — 30) 
1 1 
Dy = 107! exp [2116 (=e =) (2.563 — 0.33A + 0.026442 — 0.00067143), A>4 
1 1 
Water diffusivity in polymer phase [15] Dy = 107'° exp [2116 (= =) (—1.25A + 6.65), 3<A<4 (28) 
1 1 
Dy = 107" exp | 2416 ( — — — ] | (2.05a — 3.25), 2<A<3 
303 T 
B p 1 
Gas density in cathode p= (29) 
RT xo, /Mo, + Xwv/Mw + (1 — xo, — Xwv)/Mny 
1 
Gas density in anode pa (30) 
RT a =~ Xw)/Mp, F Xwy/Mw 
1 2.334 1 1\12 
Das = ~0,000364 (7 ) (pca pes)’ (Tea Ter)” (z + am) 
Gas diffusivity [16] p ai= v Tea Tcg Ma Mp i) 
D: = SS ———,_ D$ = Daele — s) 
RTD 4 Vi! Dix) 
€= 1 for channel e = GDL porosity in GDL, otherwise € = catalyst layer porosity 
Drag force between liquid and gas [17] MP = -MP =0.75Cp a Ug|UR (32) 
j b 
24 
Drag coefficient [18] Co = p t 0.1925 Re?) (33) 
e 
Relative velocity Ur = Ui — Ug (34) 
PgDp|UR| 
Reynold number Re = (1 — s) ———_ (35) 
Mg 
Lift force [19] My = —My = Ci pgsUr(V Ug) (36) 
Lift coefficient Cy = Cia(1 — 2.78(0.2, s)) (37) 
1/2 
Capillary pressure Pc = Pg — pi = o cosé z J(s) (38) 
p 2 3 3 4 P ys + 
Leverett TAinction tiO) i= { 1.4170. — s) — 2.12001 — s) + 1.26301 — s), ifðc < 90° hydrophilic media (39) 


1.417s — 2.120s? + 1.26383, if @¢ > 90° hydrophobic media 
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Fig. 2. Schematic view of a droplet on the GDL/channel interface. 


liquid water drop is formed on the GDL/gas channel interface. As 
mentioned above in Section 2.1, the water droplet detachment 
diameter is determined by the surface tension force and drag 
force acting on the droplet on GDL surface shown in Fig. 2. The 
surface tension force along surface direction is: 


F; = 2/0 cos cos B (=) dp (40) 


0 


Here, 04 and 6p are the advancing angle and receding angle, 
respectively, o the surface tension and 0 is the contact angle. Its 
value varies from Op to 6,4 corresponding to the angle £ from 0 to 
x. Assuming the value of 0 changes from Op to 04 periodically: 


0a — 6 

g — LOAR) (41) 
a + Or 

Integrating Eq. (40), we have 

od, g i 

Fo = -z (sin 6a + sin Op) 

1 1 
x ( ) (42) 
m —(0a— r) 1+ (Oa — OR) 


Here, de is the contact diameter between liquid and solid phase. 
For a solid surface, de can be obtained according to the droplet 
diameter d and the contact angle as: 


Oa + OR 
2: 


egies ( (43) 


But in the case of GDL in PEMFCs, the liquid water droplet is 
formed on the GDL surface by the accumulation of the water 
flowing out the GDL through the pores. So, the initial contact 
diameter between the liquid water and the solid phase equal 
the diameter of the pores. With the droplet growing, the contact 
diameter will be unchanged until coalescence happens between 
neighboring droplets. However, the liquid water droplets appear 
only at the some preferential positions but not on the whole 
surface [3], so, coalescence is not considered in the present work, 
and the contact diameter in this study is assigned to be the mean 
pore diameter in GDL. 

For a droplet on the surface in the shear flow, specifying the 
value of Op and 04 is the key problem and also is very difficult. 
They change with the gas velocity, viscosity ratio, etc. [20,21]. 
To simplify the calculation, the critical value of 6p and 64 fora 
liquid drop on an inclined PTFE plane are used for all the cases. 
For a PTFE surface [22], 


ða — 
A 299 (44) 
0A 
Assuming 6p and 04 are symmetric by the static contact angle 
Os. 
au =0.2 (45) 
bs — A0/2 ` 


For the static contact angle of liquid water on TGPH-090 sur- 
face with 10 wt.% FEP treated, Wang and co-workers [23] have 
measured it with the value of 80° under 80 °C and approximately 
100° under 20°C. 

The drag force on the droplet by the gas flow along the surface 
direction is 


F J | a a (46) 
as 2 Drez \ sing — Ow) 
m—OmM 
where 
oy = ^ > Or (47) 


Cp is the drag force coefficient, which can be estimated by 


24 


Cp = — (1 + 0.1925 Re®®?) (48) 
Re 
where Re is defined as 
d 
Ree (49) 
Mg 


And the fully developed laminar flow velocity distribution in 
a rectangular enclosure can be represented by [24] 


3 2y a 
u(y) = que i = (=) | (50) 


where v(y) is the vertical velocity, B the channel height, ve the 
average gas velocity in the channel and y represents the distance 
to the channel central line. 
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Table 4 

Boundary conditions 

Equations Inlet Outlet Wall 
Anode momentum v = 600 cm/s, u=0, p=0.2 MPa Fully developed v=0,u=0 
Anode water vapor Rey =0.7513 View =0 Vxw =O 
Cathode gas momentum ve = 700 cm/s, u=0, p=0.2 MPa Fully developed v=0,u=0 
Cathode water vapor Xwy =0.159 Viww =O Ver =0 
Cathode oxygen xo, = 0.195 Vxo, = 0 Vxo, = 0 
Cathode liquid momentum u=0; v = 0; p=0.2 MPa Fully developed v=0,u=0 
Cathode liquid saturation s=0 Vs=0 Vs=0 


In this work, gravity is neglected, so the droplet shown in 
Fig. 2 is equivalent to the droplet on the horizontal surface, and 
this is the same situation in Wang’s experiment [5]. 


3. Boundary conditions and numerical method 


Boundary conditions used in the model are given in Table 4. 
Since the model is solved by an iterative solution technique, 
the starting solution can affect the convergence and computing 
time. The starting solution for the species here was set to be 
equal to their respective inlet boundary values. The collector 
plates, membrane are impermeable for species and the GDL is 
impermeable for proton. Table 5 gives the values of parameters. 

Because the thickness of channel, diffusion layer, catalyst 
layer and membrane differ much, the non-uniform mesh grid 
is used. There are 20 mesh grids in channel, 10 mesh grids in 


Table 5 

Value of parameters 

Physical properties Value 
Faraday’s constant, F (C mol7!) 96,487 
Permeability of gas diffusion layer, Kp (cm?) 8 x 10-8 
Cathode gas viscosity, 4e (Pas) 2.0 x 1075 
Anode gas viscosity, a (Pas) 2.0 x 1075 
Liquid water viscosity, 4; (Pas) 3.565 x 1074 
Anodic transfer coefficient, œa 0.5 
Cathodic transfer coefficient, a, 0.55 
Water contact angle in diffuser, 0 (°) 120 

Gas channel width (cm) 0.05 

Gas channel length (cm) 40 

Gas channel height (cm) 0.05 
Anode GDL thickness (cm) 0.018 
Cathode GDL thickness (cm) 0.018 

Gas diffusion layer void fraction 0.7 
Catalyst layer thickness (cm) 0.001 
Catalyst layer void fraction 0.6 
Average pore size in GDL (um) 23 [25] 
Membrane thickness (Nafion® 115) 0.0127 
Cell temperature (K) 353 

Inlet pressure (MPa) 0.2 

Air inlet relative humidity (%) 100 

Fuel inlet relative humidity (%) 100 

Open circuit voltage (V) 1.17 
Operation voltage (V) 0.55 
Evap./cond. mass transfer coefficient ke, ke (kg cm~? s7!) 2 [6] 
Dissolved/vapor mass transfer coefficient, Am (s7!) 5000 [6] 
Reference kinetic parameter in cathode, Ayit, (A cm~?) 1.1 x 1076 
Reference kinetic parameter in anode, Avit (A cm~?) 1.5 x 105 
Water surface tension, o (N m7!) 0.0625 


GDL, 5 mesh grids in the catalyst layer, 20 mesh grids in the 
membrane and 450 mesh grids along the channel (y-direction). 
The governing equations were solved by our C language code 
using iterative method. The Simple Algorithm [26] is used to 
solve single phase flow in anode, and Inter-Phase Slip Algorithm 
(IPSA) is used to solve two-phase floe in cathode. When the 
relative error of all variables between two iterations are less 
than le—4, it is considered to reach a convergence. A sensitive 
analysis was conducted by doubling the number of mesh grids. 
The solution difference on average is less than 3%, so it was 
assumed to be mesh grid independent. 


4. Results and discussion 
4.1. Model validation 


To validate the model above, firstly, the predicted droplet 
detachment diameter is compared with the experimental data 
by Wang [5], which is a very important factor in the TF model 
presented in this study. Fig. 3 shows the comparison between 
calculated detachment diameter by our simplified model and 
Wang’s experimental data, and they agree pretty well. This vali- 
dates the simplified model. Second, the comparison between the 
predicted fuel cell performance and our experimental data also 
have been presented in Fig. 4 for different operation pressures. 
Good agreement is obtained, so the present model is demon- 
strated to be accurate. 
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0.030 È Calculated result for contact angle=120° 

§ 0.025 

S 

oO 

E 0.020 

& 

no} 

© 

SB 0.015 

O 

2 

a 


0.010 


0.005 


30 35 40 45 50 55 60 65 70 75 80 
Velocity,m/s 


Fig. 3. Comparison between predicted droplet diameter and experimental data. 


870 G. He et al. / Journal of Power Sources 163 (2007) 864-873 


Voltage,V 


@ 0.2MPa,experiment 

0.2MPa,calculation 
* 0.3MPa,experiment 
== Q,3MPa,calculation 


0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 


Current density,A.cm” 


Fig. 4. Comparison between numerical performance of PEMFC and experiment 
data for different operating pressures. 


4.2. Effect of GDL/channel interface hydrophobicity on 
droplet detachment diameter and water saturation 
distribution in GDL 


There are also results of predicted detachment diameter 
for different contact angle (80°, 120°) on the GDL surface 
in Fig. 3. The results show that with the increases of con- 
tact angle on GDL surface, the detachment diameter decreases. 
Figs. 5 and 6 show the water saturation in GDL for the base 
case and the contact angle on GDL surface equal to 120°. For 
both cases, the water saturation increases from inlet to outlet 
in vertical direction and from channel to catalyst layer in hori- 
zontal direction. This is consistent with the observation results 
by Yang [3]. On the GDL/channel interface, the water satu- 
ration increases from inlet to outlet. It can be demonstrated 
by the fact that with gas flows from inlet to outlet and liquid 
water moves through GDL to GDL/channel interface, more and 
more water droplets are dragged by the gas and travels with 
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Fig. 5. Saturation in GDL and channel for base case. 
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Fig. 6. Saturation in GDL and channel for contact angle on GDL/channel inter- 
face equal to 120°. 


acceleration in the channel. So the space on the GDL/channel 
interface occupied by the liquid water increases. Comparing 
Fig. 5 with Fig. 6, it is obvious that the water saturation for 
base case is larger than that of 120° case. According to Fig. 3, 
large contact angle on GDL/channel interface produces small 
detaching droplet, so small droplet travels by drag force eas- 
ily. When the droplet is small enough, the droplet will move 
at the same velocity as gas, and the efficient removal of liquid 
water can be obtained. So high contact angle on GDL/channel 
interface is advantageous for water removal in GDL and gas 
channel. 


4.3. Effect of water surface tension on droplet detachment 
diameter and water saturation distribution in GDL 


According to Eq. (40), the adhering force is proportional to 
the surface tension. Fig. 7 shows the droplet detachment diam- 
eter for normal water surface tension, 80% of the water surface 
tension and 120% of the water surface tension. The detach- 


—— Water surface tension 
80% water surface tension 
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Fig. 7. Predicted detachment diameter for different surface tension. 
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Fig. 8. Saturation in GDL and channel for 80% of water surface tension in GDL. 


ment diameter obviously increases with the increase of surface 
tension. From this point of view, reducing the surface tension 
may be a way to enhance liquid water removal on GDL/channel 
interface. However, as shown in Eq. (38), reducing surface ten- 
sion may decrease the capillary pressure in GDL, which is the 
dominating force for liquid flow in GDL. In this way, reduc- 
ing surface tension may cause liquid water accumulation in 
GDL, which is demonstrated by the numerical results shown 
in Figs. 8 and 9, where the water saturation in GDL for 80% of 
the water surface tension in GDL and 120% of the water sur- 
face tension in GDL are given, respectively, but the value of 
surface tension used to calculate the droplet detachment diam- 
eter retains original value. It is obvious that decreasing surface 
tension in GDL causes serious saturation in GDL, and this 
is disadvantageous for liquid removal. To investigate the final 
effect of surface tension (both in GDL and on GDL/channel 
interface) on the water removal, numerical simulation has been 
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Fig. 9. Saturation in GDL and channel for 120% of water surface tension in 
GDL. 
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Fig. 10. Saturation in GDL and channel for 80% of water surface tension both 
in GDL and on GDL/channel interface. 


done. Figs. 10 and 11 show the water saturation in GDL for 
80 and 120% of the water surface tension (both in GDL and 
on GDL/channel interface), respectively. It can be seen that 
water saturation in GDL for 80% of water surface tension is 
largely decreased compared with the result in Fig. 5. Along the 
GDL/channel interface, the saturation changes slightly, which 
means that the droplet moves faster for small ones, thus, only 
a few spaces on GDL/channel interface are occupied by the 
droplet. In Fig. 11, although the capillary pressure in GDL 
increases due to the increase of surface tension, the droplet on 
GDL/channel interface cannot be easily dragged by the gas due 
to the large size. It is the cause of the accumulation of liquid on 
GDL/channel interface, and correspondingly, it hinders the liq- 
uid water transfer from catalyst layer to GDL/channel interface. 
Saturation in GDL increases too. Fig. 12 shows the current den- 
sity along the y direction on cathode catalyst layer/membrane 
interface for 80 and 120% of the water surface tension. It indi- 
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Fig. 11. Saturation in GDL and channel for 120% of water surface tension both 
in GDL and on GDL/channel interface. 
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Fig. 12. Current density on cathode catalyst layer/membrane interface for 80 
and 120% of water surface tension. 


cates that there is hardly any difference in the upstream section 
for two cases, but for downstream section, due to the high 
saturation of 120% of surface tension case, the current den- 
sity is obvious lower than that of 80% case. According to the 
result above, it can be concluded that increasing surface ten- 
sion is disadvantageous for liquid removal in the GDL and 
channel. 


4.4. Oxygen mass fraction in cathode 


Figs. 13 and 14 show the oxygen fraction distribution in 
cathode channel, GDL and catalyst layer for base case and the 
operation voltage equal to 0.4 V. It indicates that the oxygen mass 
fraction decreases from inlet to outlet and from channel to cat- 
alyst layer. In addition, with the decreases of operation voltage 
(i.e. increases of current density), the consumption of oxygen 
increases, and the minimum oxygen mass fraction appears at 
the inner of downstream section. 


Fig. 13. Oxygen mass fraction in cathode channel, GDL and catalyst layer for 
base case. 
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Fig. 14. Oxygen mass fraction in cathode channel, GDL and catalyst layer for 
0.4 V of operation voltage. 


5. Conclusion 


In this paper, a two-dimensional two-fluid (TF) model is 
developed to simulate the two-phase flow behavior in PEM- 
FCs. A simplied model is also applied to predict the detachment 
diameter of the liquid droplet at the GDL/channel interface. The 
effect of the GDL/channel interface properties (contact angle, 
surface tension) on the two-phase behavior was mainly investi- 
gated, and the conclusions are: 


(1) Ahigh contact angle of the liquid water at the GDL/channel 
interface produces a small detaching liquid droplet, and 
enhances water removal at the GDL/channel interface and 
in the GDL, which is advantageous for PEMFC operation. 

(2) A low surface tension produces a small detaching liquid 
droplet, which is advantageous for water removal in the 
channel and the GDL, even though a low surface ten- 
sion decreases the capillary force in the GDL. In con- 
trast, high surface tension hinders liquid removal from 
the GDL/channel interface, which is disadvantageous for 
PEMFC operation, even though it enhances the capillary 
force in the GDL. 


In this work, the two-phase flow behavior at the anode is 
not considered, because it is not significant unless the anode is 
severely flooded. The hydrophobicity or hydrophilicity of the 
gas channel walls also have effects on water removal, especially 
for large liquid droplets, and a further study will be done in 
future. 
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